package integration;

import static mathUtils.VecMath.*;

public class Butcher extends IntegratorBase {

	float sb[][];
	float ds[][][];
	final int numTrials;  /* Number of derivatives collected INCLUDING original.
	                       * So, 1 more than what's usually called "number of
	                       * trials".  (It was either that or make the code really
	                       * ugly with a lot of "numTrials+1"'s.)      
	                       */	
	float tabA[][];
	float tabB[];
	float tabC[];
	
	public Butcher(Integrable s, int n) {
		super(s);
		numTrials = n;
		sb = vecCloneSpace(sys.state());
		ds = new float[numTrials][][];
		for (int i=0;i<numTrials;++i) {
			ds[i] = vecCloneSpace(sys.dstate());
		}
		tabA = new float[numTrials][];
		tabA[0] = null;
		for (int i = 1;i<numTrials;++i)
			tabA[i] = new float[i];
		tabB = new float[numTrials];
		tabC = new float[numTrials];
	}
	public Butcher(Integrable s) {
		this(s, 2);
		tabA[1][0] = 1f;
		tabB[0] = tabB[1] = 0.5f; // Do a Heun integrator by default
		fillC();
		printTableau();
	}

	void fillC() {
		vecSum(tabA, tabC);
		vecScaleDown(tabA, tabC);
	}
	
	public void printTableau() {
		String line = new String();
		for (int i = 0; i<numTrials; ++i) {
			line = tabC[i] + " | ";
			for (int j = 0; j<i; ++j) {
				line += (tabA[i][j] + "  ");
			}
			System.out.println(line);
		}
		line = "-----";
		for (int i = 0; i<numTrials; ++i) line += "-----";
		System.out.println(line);
		line = "    | ";
		for (int i = 0; i<numTrials; ++i) line += (tabB[i] + "  ");
		System.out.println(line);		
	}
	
	
	@Override
	public void finalizeSlope(float dt) {

		vecCopy(sys.state(), sb); // Store current state.
		vecCopy(sys.dstate(), ds[0]); // Slope 0 is from the start point, so just copy.
		
		// Run trials
		for (int i = 1; i<numTrials;++i) {
			vecZero(sys.dstate());
			vecAddScaled(sys.dstate(), ds, tabA[i], sys.dstate());
			takeStep(dt*tabC[i]);
			findSlope();
			vecCopy(sys.dstate(), ds[i]);
			vecCopy(sb, sys.state()); // Revert state after each trial.	
		}
		
		vecZero(sys.dstate());
		vecAddScaled(sys.dstate(), ds, tabB, sys.dstate()); // Construct final slope
		
	}

	@Override
	public void takeStep(float dt) {
		super.takeStep(dt);
	}
	
	public void adjustSpace() {
		sb = vecCloneSpace(sys.state());
		ds = new float[numTrials][][];
		for (int i=0;i<numTrials;++i) {
			ds[i] = vecCloneSpace(sys.dstate());
		}
	}

}
